Working with FASTQ Sequence¶
FASTQ Format Handling¶
Exploring Data
Reading Data
Cleaning Data
FASTQFormat Handling Template
# get data
!wget http://d28rh4a8wq0iu5.cloudfront.net/ads1/data/SRR835775_1.first1000.fastq
--2020-09-19 22:52:26-- http://d28rh4a8wq0iu5.cloudfront.net/ads1/data/SRR835775_1.first1000.fastq
Resolving d28rh4a8wq0iu5.cloudfront.net (d28rh4a8wq0iu5.cloudfront.net)...
13.227.214.178, 13.227.214.22, 13.227.214.185, ...
Connecting to d28rh4a8wq0iu5.cloudfront.net (d28rh4a8wq0iu5.cloudfront.net)|13.227.214.178|:80...
connected.
HTTP request sent, awaiting response...
200 OK
Length: 224786 (220K) [application/octet-stream]
Saving to: ‘SRR835775_1.first1000.fastq’
SRR835775 0%[ ] 0 --.-KB/s
SRR835775_ 14%[=> ] 31.09K 119KB/s
SRR835775_1 42%[=======> ] 92.73K 195KB/s
SRR835775_1. 96%[==================> ] 210.76K 220KB/s
SRR835775_1.first10 100%[===================>] 219.52K 229KB/s in 1.0s
2020-09-19 22:52:30 (229 KB/s) - ‘SRR835775_1.first1000.fastq’ saved [224786/224786]
Exploring Data¶
# examine first few lines
!head ../data/SRR835775_1.first1000.fastq
@SRR835775.1 1/1
TAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTCACCCTAACCCTAACCCTAACCGTATCCGTCACCCTAACCCTAAC
+
???B1ADDD8??BB+C?B+:AA883CEE8?C3@DDD3)?D2;DC?8?=BAD=@C@(.6.6=A?=?@##################################
@SRR835775.2 2/1
TAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACC
+
CCCFFFFFGHHGHJJJJJIJGIIJJJJJJJIJIJJJJJFJJFGIIIIH=CBFCF=CCEG)=>EHB2@@DEC>;;?=;(=?BBD?59?BA###########
@SRR835775.3 3/1
TAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGG
# examine first few lines
!tail ../data/SRR835775_1.first1000.fastq
+
@@CFFFFFHHHAFHI@GBHIHIIGGGE?DGGI<GGIIGHDGHIIGICHIFHHHHFDEBCECCECCCCCCC@@CACBBAB#####################
@SRR835775.999 999/1
GGCTGCCCCTGCAGGTCCCTGTCACCTCTCACATGTCCCTGCCTAATCTTGCAGGTCCCAGAGAACTACTTCTATGTGCCAGACCTGGGCCAGGTGCCTG
+
@@CFFDFFHHFGFJJGHIICHEHGJJIEGIIJJJJJJJJJJJJJIGGJJJIJJJIHHGHIJJJIHIIJJHGHHHGHFFFFFFCDDDCBDDDDDD@CCCDC
@SRR835775.1000 1000/1
GTTGGCTGCCCCTGCAGGTCCCTGTCACCTCTCACATGTCCCTGCCTAATCTTGCAGGTCCCAGAGAACTACTTCTATGTGCCAGACCTGGGCCAGGTGC
+
C@CFFFFFHHHHGHGGHIFIIJJJHJJJIGIIJIIIIIIJIJJJJJEBGIJJJJJGGJBFEFGHIJEGEEEHGECBDFFFFEEEDDDDDC=BBDDBB>AD
# convert quality score character to numeric
ord("?")
63
ord("B")
66
ord("#")
35
Reading Data¶
# open file
with open("../data/SRR835775_1.first1000.fastq", "r") as f:
seq = f.read()
seq[:500]
'@SRR835775.1 1/1\nTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTCACCCTAACCCTAACCCTAACCGTATCCGTCACCCTAACCCTAAC\n+\n???B1ADDD8??BB+C?B+:AA883CEE8?C3@DDD3)?D2;DC?8?=BAD=@C@(.6.6=A?=?@##################################\n@SRR835775.2 2/1\nTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACC\n+\nCCCFFFFFGHHGHJJJJJIJGIIJJJJJJJIJIJJJJJFJJFGIIIIH=CBFCF=CCEG)=>EHB2@@DEC>;;?=;(=?BBD?59?BA###########\n@SRR835775.3 3/1\nTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGG'
# remove name line
with open("../data/SRR835775_1.first1000.fastq", "r") as f:
seq = f.readline()
seq = f.read()
seq[:500]
'TAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTCACCCTAACCCTAACCCTAACCGTATCCGTCACCCTAACCCTAAC\n+\n???B1ADDD8??BB+C?B+:AA883CEE8?C3@DDD3)?D2;DC?8?=BAD=@C@(.6.6=A?=?@##################################\n@SRR835775.2 2/1\nTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACC\n+\nCCCFFFFFGHHGHJJJJJIJGIIJJJJJJJIJIJJJJJFJJFGIIIIH=CBFCF=CCEG)=>EHB2@@DEC>;;?=;(=?BBD?59?BA###########\n@SRR835775.3 3/1\nTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGG'
# read base sequence
with open("../data/SRR835775_1.first1000.fastq", "r") as f:
seq = f.readline()
seq = f.read()
seq[:500]
'TAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTCACCCTAACCCTAACCCTAACCGTATCCGTCACCCTAACCCTAAC\n+\n???B1ADDD8??BB+C?B+:AA883CEE8?C3@DDD3)?D2;DC?8?=BAD=@C@(.6.6=A?=?@##################################\n@SRR835775.2 2/1\nTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACC\n+\nCCCFFFFFGHHGHJJJJJIJGIIJJJJJJJIJIJJJJJFJJFGIIIIH=CBFCF=CCEG)=>EHB2@@DEC>;;?=;(=?BBD?59?BA###########\n@SRR835775.3 3/1\nTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGG'
seq = seq.replace("\n", "")
seq[:500]
'TAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTCACCCTAACCCTAACCCTAACCGTATCCGTCACCCTAACCCTAAC+???B1ADDD8??BB+C?B+:AA883CEE8?C3@DDD3)?D2;DC?8?=BAD=@C@(.6.6=A?=?@##################################@SRR835775.2 2/1TAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACC+CCCFFFFFGHHGHJJJJJIJGIIJJJJJJJIJIJJJJJFJJFGIIIIH=CBFCF=CCEG)=>EHB2@@DEC>;;?=;(=?BBD?59?BA###########@SRR835775.3 3/1TAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGT'
seq = seq.replace("\t", "")
seq[:500]
'TAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTCACCCTAACCCTAACCCTAACCGTATCCGTCACCCTAACCCTAAC+???B1ADDD8??BB+C?B+:AA883CEE8?C3@DDD3)?D2;DC?8?=BAD=@C@(.6.6=A?=?@##################################@SRR835775.2 2/1TAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACC+CCCFFFFFGHHGHJJJJJIJGIIJJJJJJJIJIJJJJJFJJFGIIIIH=CBFCF=CCEG)=>EHB2@@DEC>;;?=;(=?BBD?59?BA###########@SRR835775.3 3/1TAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGT'
FASTQ Format Handling Template¶
def readFastq(filename):
"""Reads FASTQ file and remove the special characters!"""
sequences = []
qualities = []
with open(filename) as fh:
while True:
fh.readline() # skip name line
seq = fh.readline().rstrip() # read base sequence
fh.readline() # skip placeholder line
qual = fh.readline().rstrip() #base quality line
if len(seq) == 0:
break
sequences.append(seq)
qualities.append(qual)
return sequences, qualities
seqs, quals = readFastq('../data/SRR835775_1.first1000.fastq')
seqs[:20]
['TAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTCACCCTAACCCTAACCCTAACCGTATCCGTCACCCTAACCCTAAC',
'TAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACC',
'TAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGG',
'TAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAAGGGTTGGGGGTTAGGGGTAGGGGTAGGGTTA',
'CTAACCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCAACCCTAACCCTAACCCTCACC',
'AACCCTAACCCTAACCCTAACCCTAACCCTAACCCCTAACCCTAACCCTACCCCTAACCCCCAACCCTCACACCAACCCTAACCCTACCCCCAACCCCAC',
'TAGGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGCTGGGTTAGGGGTAGGGTTAGGGTTAGGGTTAGGGGTAGGAGTTCGGGAGAGCACACG',
'TAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCAACCCCAACCCCAACCCCAACCCTAACCCTAAACCCAAACCTAAA',
'AGGGTTAGGGTTAGGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGGTTAGGGTTGGGGTTAGGGTTGG',
'TAGGGTTAGGGTTAGGGTTAGGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTGGGGTAGGGGTTGGGGTAGGGGTAGGGGTGGGGGTGGG',
'TAACCCTAACCCTAACCCTAACCCTAACCCAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTACCCTACCCTAACCCT',
'AGGGTTAGGGTTAGGGTAGGGTTAGGGTTAGGGTTAGGGGTAGGGGTGGGGGTAGGGGGTAGGGGGTGGGGGGAGGGTTAGGGTTTGGGTTCGGGGTAGG',
'ACCCTAACCCTAACCCTAACCCTAACCCAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCTACCCA',
'CCTAACCCTAACCCTAACCCTAACCCAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTACCCTACCCTAACCCTAACC',
'CTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAA',
'AACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACC',
'AACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCTTAACCTAAACCCTAACCCAAACC',
'GGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTAGGGGTAGGGGTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGGTAGGGTTA',
'GTTAGGGTTAGGTTAAGGTTAAGGTTAGGGTTAGGTTAAGGTTAAGGTTAGGTTAAGGTTAAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGG',
'TTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTAGGGTTCGGGTTAGGGTTAGGGTTAGGGTTAGGGTTGGGGTTGG']
len(seqs)
1000
quals[:20]
['???B1ADDD8??BB+C?B+:AA883CEE8?C3@DDD3)?D2;DC?8?=BAD=@C@(.6.6=A?=?@##################################',
'CCCFFFFFGHHGHJJJJJIJGIIJJJJJJJIJIJJJJJFJJFGIIIIH=CBFCF=CCEG)=>EHB2@@DEC>;;?=;(=?BBD?59?BA###########',
'@@<DD?DDHHH<CBHII:CFGIGAGHIIG?CCGGE0BDHIIHIGICH8=FHGI=@EHGGGEEHH>);?CC@.;(=;?59,5<A599?CB>ABBCB><88A',
'@CCFFDDFHHHDFHIJJCGGIJJHIIHJC?DHHIJ0?DGHI9BBFHICGGIGI=CDEGI=?AAEF7@?################################',
'@<@FDFDDBBFHH@EBGGGEH@FFHE;EHIEGGBHG?)9?8BDGGBGGBCDGI=93=C6==C;CCD(?@>@#############################',
'CCCFFDFDHGHHDHHGIIIHGIJ<GBBBFBGGHDIH@B)?)0?F>DB#####################################################',
"???=DD0@3ADB<<CE;+3<;EE?1:)1?C1?3BDCB@;)99)0;@;)-(.7A'55A1?.6;(6;.5;;??#############################",
"@@@DDF?EFFHA;DDHGI<FDFGECCHF:C9DH9F3):??GHB)BFHEAHBG;C';C@EE45?#####################################",
'@@@FADDFHGBFFHJJJFGEHGGGGCEGHCBFGGJ:D?FGIDHDFGIBFDGIJ=CEHDH7?=BFD@A;>AB5=2=?BB######################',
'@CCFFADDHHH<CBGIICGGHIJI?FFGHJ?DHGIJ?DGIIIBFCHHI8@@GEE.==A?C########################################',
'CCCFFFFDHDHHHCIEIJJJICEGHIJIJIHIIII>FFGGG>GCGH9(B@FGE>EEGCH>A;C(;;;3>>A2;(;;A99(<A<(9ABD@?##########',
'<@<BDDDDD8CD?CGE<<CA;<<+2:C9:??FC*:8BGG)??##########################################################',
';?@FFDADFHGBHIGGHIGGGGGI;FGCD8CFH3B?BD89?)?FH>;=C;C;==@DHGH<EEED;)6>AB@(=A(5;99?B###################',
'CCCFFFFFHHGHHJIIJIJJIIIJJJFHHHHGCHGIJGIJGDHEIJIGIDHH@CH;DEEFC6?BB@D(;>;=A(;5=99(59?B################',
'C@CFFFFFHHHHHIIJJJJGJIIJJIIJEHGGHJGIIEI@FHIIIJIFCGHCEHCEECH>E@CEDDFBECCCDDCCB?A?ACDBD<<3<8<AC?B8?8<9',
'<<@BDDAB?CFDH;CCE;BAGHICHGHEHIIIIIG;?FFGHIAHIHHGIIIIIB;=@EDD>AEH=D>>;@;@>===A(=CBBB@??99?<<?8?2?<C?B',
'CCCFFFFFHHHHHIJJJJJJJJJJJJJHGGIIJJJJJJJJJJIIJIJJJJJGGIJJIIJJEHEHEFEEFE2;@ABD(;@CDD##################',
"BBBDDFFFDHHHHJEHHIJJCGEHIJ:EFIHI?GDHHI9BFHIIB?;CF'B;FHI'=:ADF.;;>>@(,(55?,559?A(5(8<B(8289?&0)2??322",
'@@@FFFDDFDHFFHEIIGHG@EHJFGGGGCEGCGCCFBEGI?FHGIJDDFGH?FC@DH8B@:CAGCAEGIG@E2A?D7;BEDE;@=?C?;9?BD?5<?BC',
'B@CFFFBDHHHHAGGFHJCGEHIJFHGIHI?HFGIJ:DDGHJDHFFHI=FCHHG@EAEHHAE9@DD?B;ACD;?(9<A>A9<AB+82<CD##########']
ord('#')
35
ord('?')
63